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o : 

Analysis of displacement in a Hele-Shaw cell and porous media is a source of a multitude of 
mathematical problems which provide some insight into general features of nonlinear boundary 
dynamics([4],[8],[l]). Here, we consider a slightly modified version of the classical problem of flow 

7— I ■ in a potential external field which displays some new features related to existence of singularities 

of the external field. The study was prompted by the interest in coupled flow phenomena in 

y— I | saturated porous media in the presence of electric current ( "electrokinetic phenomena" ) . This 

specific case will be briefly discussed later. However, eventually it became clear that the topic 
■ deserves study per se. Throughout the text we will use expressions 'Hele-Shaw flow' and 'flow 

through a porous medium' interchangeably as synonyms due to the well known analogy between 
the Darcy law for a porous medium and the flow rule in a thin gap between two parallel walls. 

m : 
o ■ 

^ : 1 The Problem Statement 

We assume that fluid occupies a finite domain D(t) in the plane (x, y) of a Hele-Shaw cell 
of gap thickness b. The flow velocity w is determined by the flow rule 

" . k„ k 

>■ w = Vp+-g. 1.1 

Here, p is the fluid pressure, p is the fluid viscosity, k = b 2 /12 is the gap "permeability", 
g = {g x ,g y } is the body force due to an external field. We assume the body force to be 
potential, 

g=-W(z,y). (1.2) 
Two familiar examples are the gravity force and the centrifugal force, corresponding to 

* = pgh, and * = --pu 2 r 2 , (1.3) 

respectively. Here, p is the fluid density, h is the height above a datum, g is the acceleration 
due to gravity, r is the distance from the rotation axis normal to the cell plane; u is the 
rotation rate. These cases allow thorough study; they are considered in particular in [4], [3]. 
Here, we are going to study the case when the external potential \I/ has singularities within 
and/or outside the domain D(t). 



1 



The flow field satisfies the continuity equation 

N 

V ■w = ^q j 5(x-x j ,y-y j ). (1.4) 
j'=i 

Here, are strengths (flow rates) of the point sources (sinks) within the flow domain. 
We assume that at the boundary T = dD(t) the pressure vanishes, 

p(x,y) = Q, (x,y)eT. (1.5) 

The boundary dynamics is governed by the relation 

v n = w n = w-n, (1.6) 

n being the outward normal to T, and v n velocity of propagation of the boundary in the 
normal direction. 

We assume now that the external potential field satisfies the equation 

M h 
AG = Y,Q m 8(x-xf m ,y-y' m ), G = — tt, (1.7) 

m=l ^ 

in the entire plane (x, y) with boundary condition 

|VG|->0, (x 2 + y 2 )^oo. (1.8) 



Therefore, 



M 

G — ReF(z), F = Y,^Hz-*L), 

m=l 



z = x + iy; z' m = x' m + iy' m . (1.9) 
Let us introduce the velocity potential 

$(x,y) = --(p + V). (1.10) 

It satisfies the following problem: 

TV 

A<f> = -^2q j 5(x-x j ,y-y j ), ( Xj , yj ) G D(t); (1.11) 

$(x,y) = --*(x,y) = G(x,y), (x,y) G T(t); (1.12) 

Wn = — , (x,y)er(*). (1.13) 

The last equation serves to describe the moving boundary dynamics. 
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The only difference with the usual Hele-Shaw problem is that the flow potential does 
not vanish at the boundary, but should be equal to a specified function of the boundary 
point. 

Let now u(x, y) be a harmonic function in a domain D* such, that D(t) remains within 
D*. Then as a straightforward generalization of the Richardson Theorem [10], [11], [4] we 
find: 

d f N f 

— udA = y2lM z j) + VG-VudA. (1.14) 

dt J D (t) j=1 JD(t) 

The following chain of equalities proves statement (1.14): 

d f f <9$ „ 

u-—dl = 



[ udA= [ 

Jdu) Jdi 



dt jD(t) JdD(t) on 

du f 

dl + / (wA$ - $Aw) dA 



8D(t) 9n j D ( t ) 

r du N r N 

/ G 7^ dl + E = V • {GVu)dA + 

JdD(t) on ^ J D{t) ^ 

N 

i=l ^( 



The l.h.s. of Eq.(1.14) is a time derivative of a harmonic moment of the domain D(t). 

This equation leads to explicit analytic techniques of predicting domain evolution 
provided that the operator VG • V maps harmonic functions to harmonic ones (and the 
domain initially belongs to a certain class of domains). It can be shown, that it is possible 
only, if 

G = ax + by + c(x 2 + y 2 ) + d, (1.15) 

with constant a, b, c, d. Essentially, it is a combination of uniform ("gravity") and axisym- 
metrical ("centrifugal force") fields treated previously [3], [4]. Of course, for c ^ 0, it can 
be reduced to pure rotation about a shifted axis. 

Here, we will be interested primarily in equilibrium shapes of the flow domain under 
combined action of the flow and the external potential. In this case, one can derive 
effective solutions for a much wider class of the external potentials. 

2 Steady-state shapes 

For the steady state (equilibrium) domain D we arrive at the following moment problem: 

Vm : Am = in D, (2.1) 



L 



N 



VG ■ VudA = -J2 a M z j)- (2-2) 
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(Of course, equilibrium domain can exist only if the net fluid flux vanishes, X^jli Qj = 0-) 
Now we assume that the external field has the form 

M 

G(z) = J2^\z-z' m \- (2-3) 



m=l 



In other words, it can be considered as an electric field generated by a finite array of point 
charges in a plane. For brevity sake, we will refer to it as an "electric potential" field. 
Note that z' m can be both inside and outside D. Let us now introduce the corresponding 
complex potential F(z) and 'complex current' uj{z): 

M 

F{z) = G(z) + = X < z - ( 2 - 4 ) 

m=l 

W ( Z ) = = ^ - ^ = V - Qm - . (2.5) 

(The complex potential -F(z) is, generally speaking, multivalued, unless all the 'electric 
sources' are outside D.) Then the moment equation (2.2) can be written as 

r N 

J D = / Lu(z)U'(z)dA = -J2^U( Zl ), (2.6) 
Jd i= i 

for an analytic function [/ on a neighbourhood of D. 

Integral J D in the l.h.s. of Eq.(2.6) converges even if some z' m belong to D as oo(z) has 
simple poles at these points. We write the integral as 

J D =fZ^jU'(z)dxdy=^(fZ^jU(z)dz+ %tfOO- (2-7) 

The last transformation follows from the Green Theorem; the sum in the r.h.s. accounts 
for contributions of poles of u in D. 
Let us now choose 

U(z) = 1 , 

71{W — Z) 

w being a point outside D. Then Eq.(2.6) becomes 

2mJ 9D w-z ^n(w- Zi ) m .^ eD 2?r w - z' m 
We denote the primitive of the function in the r.h.s. of Eq.(2.8) by h(w): 

M«0=£fln(«;-*)- E ^Hw-z'J- (2.9) 

i=i m: z^er* 
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The following Theorem due to Richardson [10] plays a crucial role in solving the problem 
of finding the domain D: 

Theorem 

Let / : K — > D be a conformal mapping that maps unit disk of the £-plane onto D. 
Then the function 

4? (nnb) - Hno)) (2.10) 



initially defined in a vicinity of the unit circle extends analytically to a holomorphic 
function in K. 

Proof. On the unit circle ( = e lip , ( = l/( and thus it suffices to show that the 
differential d(F(z) — h(z)) on dD extends to a holomorphic differential in D. 
So according to the Cauchy Theorem it is necessary to check that 



= for tiD . (2 . n) 

dD t — Z 

However, it follows directly from Eqs.(2.8),(2.9) and the fact that dh is holomorphic 
outside D. 

This implies the following important corollary: 
Corollary. 

The function J^F(/(£)) is rational. 
Proof. For any function f(z) denote 



/*(*) = m- 

It implies immediately 

WW) = F*(r(z))- 

Then according to the Theorem (2.10) the differential 

d(F'(/'(i/c))-M/(C))) 

is analytic in the unit disc. 
Therefore, the differential 

D(C) :=d(F*(r{±)) 

has the same singularities as dh(f(()) in the unit disk K, i.e. simple poles at the points 

c = rVm), and c = r 1 fe), j = l,...,N. 

Therefore dF(f(Q) has poles at the points 



l/f-^zj), l/f-^z'J, z' m eD 
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outside the unit disk. On the other hand, dF(f(()) has poles at the points 

C = r 1 (4J, z' m eD 

within the unit disk. 

Therefore, S)(C) has a finite number of poles, and hence it is rational. Then determi- 
nation of the precise form of f(() can be reduced to a set of nonlinear algebraic equations. 

Note, that this analysis can be in a standard way extended on the limiting case of 
coalescence of hydrodynamic sources and sinks corresponding to multipoles. In such a 
case, the terms 

ir(z - zj) 

have to be replaced with the terms 

, n > 1, etc. 



3 Harmonic potential 

3.1 Univalent F(z) 

Assume that F\z) ^ 0, oo, and that F(z) is univalent in D. Then we can solve the 
problem in a more straightforward way. We just notice that for the domain D = F(D) 
one has 

[U'(z)dxdy = -J2^U(F 3 ); ~z = F(z); F 3 = F( Zj ). (3.1) 

JD 

However, this is exactly the form of the moment equation that corresponds to uniform 
external field with G(z) = x (the Hele-Shaw flow in the presence of gravity) [4]. If we 
take u(z) = z n , Eq.(3.1) becomes 

M n _ x = / ~z n ~ l dxdy = - V qjFf/n, (3.2) 
Jd 

so that all moments are specified at given F(z). 
Letting n = 1 we have 

S = -Y,q j F j >0. (3.3) 

It is a necessary condition for the existence of a steady-state solution. In particular, 
if the flow is generated by a dipole, (a source-sink doublet of strength ±q = ±/i/(2e) at 
z = ±e) then as e — > the r.h.s. of Eq.(3.3) tends to /xF'(O), and Eqs.(3.2) becomes 



= / dxdy = n\ M„ = 0, n = 1,2,3, ... (3.4) 
Jd 
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Obviously, in the F-plane the equilibrium domain is a circle of the radius 



Thus, there is a fixed value of the area of the equilibrium domain in the potential plane 
for which such a domain exists. Obviously, the equilibrium in this case is due to a fine 
balance between hydrodynamic and external forces, and is unstable. 

The above elementary example is generic in the sense that for given set of hydrody- 
namic and electric sources the area of the equilibrium domain, provided such a domain 
exists, can assume only a discrete set of values, if the electric sources are outside the 
domain, so that F(z) is analytic in D. 

It is reasonable to ask about the fate of a domain that evolves under combined action of 
the balanced hydrodynamic sources and the external field starting from a non-equilibrium 
shape. While in general case the answer is beyond our capacities, some insight can be 
derived from the simple case of "gravity" , i.e. uniform potential field, 

F(z) = gpz. (3.5) 

Then the moments dynamics equation (1.14) becomes 

4- [ udA = Y q lU { Zl ) -C [ ^dA; C = ^. (3.6) 

If we take now 

1 „, „ . . f dA 



u(z) = —^ -, weZ\D; xW= / 

n(w - z) J D 

then Eq.(3.6) can be written as 



tt(w — z) ' 



dx(z,t) dx(z,t) ft 

dt dx ^^(z-Zi)' { ' 

It is a first order p.d.e. that is readily solved explicitly. The solution has the form 

r* N ndr 

«'•<> = Uz + ct) + 1 g * + cl-T)-*y (3 - 8) 

For the dipole, Xo( z ) — A/z, and the integrand becomes fi/[n(z+C(t — r)) 2 ], and therefore 

u AC — a 



cz c(z + cty 

The first term is the Cauchy transform for a circle of area Aq = fi/C centered at the 
origin; the second term corresponds to the circle of the area A — A centered at z — —Ct. 
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It means that at large t the solution is combination of a steady-state circle of area A at 
the origin corresponding to the limiting steady-state solution, and the circle of the area 
A — A "sinking" in the gravity field. 

As general solution (3.8) is valid for any function h, it is tempting to state, that an 
arbitrary initial domain of sufficiently large area under the combined action of gravity and 
a dipole at the origin (z = 0) eventually splits into two parts, namely, a stationary disk 
of area A = fi/C centered at z — 0, and a "sinking" domain with the Cauchy transform 
of the form 

X i = Xo(z + tC)-fi/[C(z + tC)]. (3.9) 
At t — the shape of the "sinking" domain is specified by the Cauchy transform 

h 1 (z,0)=h (z)- f jL/(Cz). 

It can be derived from the initial domain in the absence of gravity by placing a sink at the 
origin and sucking the amount of fluid corresponding to the equilibrium domain area Aq. 
If we now allow this new domain to slide far enough in the gravity field, and then inject 
back the same amount of fluid at the origin, we will get exactly the Cauchy transform 
specified by the Eq.(3.9). Now we see, that this conjectured form of evolution will indeed 
occur, if the initial domain will evolve smoothly during the initial sucking of fluid. It 
will be certainly so, if the initial domain itself can be produced from a simply connected 
domain by injecting the amount of fluid Aq without violating the simply-connectedness 
condition. 

This argument allows one to develop a number of explicit solutions for domains evolv- 
ing in gravity field in the presence of a dipole at the origin. 

Now we are going to consider some less trivial examples of equilibrium domains. It is 
worth noting that such domains can exist only for special combinations of hydrodynamic 
sources and the external potential. Indeed, if we let in the moment equation (2.6) U = 
F(z), it results in 




and since the l.h.s. of this equation is positive, the (complex) electric potential should 
satisfy the condition 

N 

"I>n*;)>0. (3.10) 

A priori the sum in the r.h.s. can be any complex number, but for the equilibrium 
shape to exist, the r.h.s. should be real and positive. Obviously, this inequality can hold 
only for special form of potential. Say, in the case of a doublet source-sink of equal strength 
both of them should lie on the same force line of the electric field (ImF(zi) = lmF(z 2 )). 

As we will see later, if some electric field sources ("charges") are within the flow 
domain D, there is a continuous spectrum of the equilibrium domain areas. For example, 
in the simplest case of absence of hydrodynamic sources, qj — 0, there is no flow within 
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the equilibrium domain, the potential $ = const in D, and boundary condition (1.12) 
implies that G = const along dD, so the boundary should be a level curve of the electric 
potential. Say, in the case of a single electric source the equilibrium domains are circles of 
arbitrary radius centered at the source. Notice, that since the boundary of the equilibrium 
domain in the absence of hydrodynamic sources is just a level curve of electric potential 
the analytic continuation procedure reduces to the reflection principle of electrostatics 
applied in the C-plane. 

Example 1. Let we have two hydrodynamic sources of the strengths q\ = —q<i — q at 
Z\ = a > 0, z 2 = b respectively, and an electric "charge" Q at z = 0. Then 

F(z) = ^\nz, (3.11) 

Inequality (3.10) implies that qQ\n(b/a) > 0, so that b > a for positive Q. 

Using reduction to the "gravity" case technique, we have in the potential plane z = 
(Q/2tt) Inz the Cauchy transform for the transformed domain D = F(D) 



q 



z — In a 



Then the conformal map / of the unit disk K on D is given by the expression (cf. [3]) 

/(C) = -m^4 + f InVab, (3.13) 
with a determined from the equation 



«1b1±4 = ^1baA 0<a<l. (3.14) 
7T I -a 2 2ir Va 



Therefore, 



and 



' M (6/ a )A/2 + i' 2q' 



/(C) = e% m = Va~b [y^) ■ (3-15) 

The solution is a circle for A = 1, i.e. Q = 2q. The solution remains physically sensible 
only for not too large b/a; otherwise overlapping of different parts of predicted domain 
occurs. 

Consider condition of non-overlapping of the mapping (3.15). The condition of overlapping 

is 

/(0 = /(0 = /(VO; C = e i4> 



l/A 
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In our case it implies 



1-aCJ \l-a( 



--1 



In the plane 

_ l + a( 

U — i 7 

1 — aC 

the boundary is a circle of the radius n = 2q/(1 — a 2 ) centered at (1 + a 2 )/(l — a 2 ). As 

z = C(l /X , 

the overlapping occurs at the point of maximum arg(£i). But 

(3 = max(arg(Ci)) = sin _1 (2a/(l + a 2 )). 



So overlapping occurs at 



or 



f3 = ir\; 2a/(l + a 2 ) = sinvrA 



ttQ (b/a) x / 2 - 1 2 vrQ 

a = tan — — ; — —7- = tan — — 

4q (b/a) x / 2 + 1 4g 



(a/6) = cos ^ . (3.16) 



2gy 

Therefore, non-overlapping equilibrium domain exists in the range of parameters 

1 > a/b > (cos^A) (2/A) ; A = Q/(2q). 

For small A non-overlapping equilibrium domains exist only in a narrow range of b/a close 
to unity. If we let the ratio b/a tend to the critical value (cos7rA) 2 / A , the area of the equilibrium 
domain increases rapidly, and the domain acquires horseshoe shape. For A = 1/2 the equilibrium 
domain remains simply connected for any a/b. In the limiting case 

/(c) " (^?) 2 

the critical domain is the entire z plane with a cut along the negative real axis. 

An example is presented in Fig.l. For curve 1 the ratio Q/2q is quite close to the 
critical value. Consider now a limiting case when the hydrodynamic source and sink 
collide, and form a dipole of the moment /i. Formally, it corresponds to 

b/a = e € ; q — —, e — > 0. 
ea 



Then Eq.(3.15) becomes 
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Figure 1: Equilibrium domains. Flow is driven by a source at z = a, and sink at z = 
b, an electric charge Q is located at z = 0; plots correspond to q=l; a=l; b=4; and 
Q = 0.2734,0.2959,0.3189,0.3424,0.3664,0.3909. for curves 1-6 respectively (a). Figure 
(b) shows blow-up of the upper figure illustrating that the electric charge is outside the 
equilibrium domain. 
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So we find 

/„(C) = Jim MC) = hm „ (l + ^ J "° = « exp ( C ^f) . (3.19) 

This mapping corresponds to a non-overlapping domain iff 

\Z2\ijaQ < tt; fi/aQ < n 2 /2. (3.20) 

At greater values of fi/aQ, there is no simply connected equilibrium domain. It can be 
conjectured that in this case the electric field is too weak to prevent breakthrough caused 
by the hydrodynamic dipole. 

Figure 2 shows shapes of equilibrium domain for (j, — 1, a — 1; Q— 0.2026; 0.2410; 
0.2866; 0.3408; 0.4053. 



3.2 Singular points technique. 

Example 2. Consider now the equilibrium domains corresponding to a hydrodynamic 
dipole of the moment /i and an electric source of the strength Q both located at z — 0. 
In this case, the differential d[F(f(Q)] has singularities only at and oo. The general 
method described above implies that 

d[F(f(Q)]=(j + R)dC, (3.21) 

F(f(()) = P\n( + R( + C, (3.22) 
/(C) = Ae 2 ^ Q = At 2 * p l Q e^ R l Q K. (3.23) 
Since f(() is a conformal map, /'(0) ^ 0, oo. Therefore, 

P=^; f(() = ACe B <. (3.24) 
In order to find B, we use the the moment relation 

/ u{z)U'{z)dS = fiU'(0). (3.25) 
Jd 



In our case it assumes the form 



L 



Q -U'(z)dS = fiU'(0). (3.26) 



2ttz 
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Figure 2: Equilibrium domains. Flow is driven by a dipole of moment \x at z = a, an 
electric charge Q is located at z = 0; plots correspond to [1=1; a=l; and Q= 0.2026; 
0.2410; 0.2866; 0.3408; 0.4053. for curves 1-5 respectively (a). Figure (b) shows blow-up 
of the upper figure illustrating that the electric charge is outside the equilibrium domain. 
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Choosing U = z, we get 

f d S 27r/i . . 

J B T=-i- (3 27) 



Substituting 

we get 



z = A(e B <, 
\A + AB(\ 2 \e B (\ 2 



da 



nA + ABQ(A + ABQe« ^ 
J K AC Q 1 ; 

Evaluating the integral in the l.h.s. of Eq.(3.28) we find 

B = ^, f(()=A(e&. (3.29) 

The parameter A still indeterminate is a size parameter. It should satisfy the condition 
that /(C) is single-valued. Therefore critical values of the parameter correspond to 

(1) /'(C) = 0,_ C G dK, or 

(2) /(C) = /(C) = /(C 1 ), C G dK. 
The first condition results in: 

/'(C) = A(l + B()e B < = 0, C G dK, =► 5 = 1, or || = 1. (3.30) 
The second condition implies 

1 = /(C) / = ^^(C-C- 1 ) = e 2^+2^ sin0 (3 31) 

Then the critical condition becomes 

2/i 

L(0) = + sin = vrA;, < < vr. (3.32) 
As L(7r) = 7r, a root in the segment (0, ir) appears as derivative 

becomes negative, or at 

2 ^ 1. (3.33) 



AQ 

This condition coincides with (3.30). Hence a simply connected equilibrium domain exists 
for 

% <1. ^>|. (3.34) 
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Figure 3: Equilibrium domains. Flow is driven by a dipole of moment ji and an electric 
charge Q located at z = 0; plots correspond to /i=l; Q=l; and A = 1, v2, 2, 2 a/2, and 4 
for curves 1-5 respectively. Self-intersecting boundaries 1 and 2 correspond to non-physical 
domains. 



Thus there exists a continuous spectrum of sizes of equilibrium domains that is bounded 
from below. 

This result can be interpreted as inability of a given "charge" to prevent breakup due 
to action of hydro dynamic dipole, if the domain area is too small, or the domain boundary 
is too close to the dipole. 

Figure 3 shows boundaries of the equilibrium domains desribed mapping of the unit 
disk given by Eq.(3.29) for Q = 1, /i = 1, A = 1,^,2,2^, and 4 for curves i-o 
respectively. Example 3. Now we consider interaction of a hydrodynamic quadrupole at 
the origin with two "charges" of strength Q at z = ±a outside the equilibrium domain 
D. We fix the conformal mapping of K on D by conditions /(0) = 0; /'(0) > 0. In this 
case F(f(Q) has a pole of second order at infinity, and hence 




(3.35) 



and 




(3.36) 
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Therefore 



/(£) = z = a\jl - e ~^ e , a>0. 



(3.37) 



The parameter a depends on the strength (3 of the hydrodynamic quadrupole, namely, 
the pole coefficient of F(f(()) at oo is the same as that of h(f(Q) at £ — 0. As 



h = 



P 



the requirement implies 



P 



2na 2 ^ 



= a: a 



2ltZ 2 ' 

2_ PQ 



4n 2 a 2 ' 



Therefore the solution exists for (5 > 0, and 



a 



27ia ' 



Thus, 



/(C) = a 



\ 



1 — exp 



P 2 



a 2 Q 



C 2 = < 



\ 



This mapping remains univalent until 



exp 



a 2 Q 



remains single-valued, so that 



(3.38) 



(3.39) 



1 P , 

ay Q a 2 Q 



< 7T . 



Figure 4 shows mapping of the unit disk given by Eq.(3.39) for (3 
0.2251,0.2677,0.3183 = 1/tt, 0.3785 (curves 1-4 respectively). 



(3.40) 



1,Q = l,a 



4 Non-harmonic External Field. Reduction to the 
Riemann-Hilbert problem 

The problem of finding stationary shapes of flow domain in an external field can be ap- 
proached in a different way that allows extension to non-harmonic external field potential 
of special form. 

Let D be an equilibrium domain for a specified set of hydrodynamic sources and 
multipoles corresponding to logarithmic singularities and poles of the velocity potential 
W(z) in the field of the external force potential G(x,y). 
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Figure 4: Equilibrium domains. Flow is driven by a quadrupole of moment (3 and 
two electric charges Q located at z = ±a; plots correspond to (3—1; Q=l; and 
A = 0.2677; 0.3183 = 1/tt; 0.3785; 0.4502 for curves 1-4 respectively (a). Blowup 
of Fig.4,a. The charge remains outside the equilibrium domain; only non-intersecting 
boundaries correspond to physically admissible equilibrium domains. 



Then W(z) is an analytic (while may be multivalued) function in D having prescribed 
set of singularities; its differential dW(z) is a meromorphic function in D and 

W(z)\ zedD = G(z,z). (4.1) 

Here, 

G(M)eG^^)1( z -z)|. (4.2) 

We are going to show that in number of cases due to special form of the potential 
G(x,y) it proves to be possible to find the domain D explicitly. Consider once more the 
conformal map / : K — > D and define 

e(0 = w(/(0), (eK. (4.3) 

This function is analytic up to singularities of the specified type (poles and logarithmic 
singular points) in K and assumes real values along the boundary OK. Then it can be 
analytically continued into the entire complex plane ( using the symmetry principle: 



9(C) = ©(=). ICI>1- (4-4) 



Then 

N 



(o = E f^ ln (c - o) + in(c - c; 1 )]- (4.5) 

Hence 6(C) is known in the entire complex plane ( up to locations of the singularities Q. 
At the boundary of the unit disk 

&(() = G(z,z). (4.6) 
For given Q, it is an equation for the conformal mapping f(() that can be written as 

G(/(0,/*(i/0) = e(c). (4.7) 

In general, it is not clear how to determine the conformal map f(() from this equation. 
However, it proves to be possible under some additional assumptions on G. 
Some of these particular cases are presented below. 

4.1 Harmonic velocity potential 

Suppose that G is a harmonic function with maybe a finite set of logarithmic singular 
points within D, and let F(z) be respective complex potential. Then 



G(z,z) = -(F(z) + F(z)) = -(F(z) + F*(z)). (4.8) 
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Introducing this expression into Eq.(4.7), we get 



*"(/*(V0) = 20(C) -F(f(()). (4.9) 

Therefore F*(f*(l/()) has finitely many singularities within K. Those inside K cor- 
respond to singularities of 6(C) and F(f(()) inside K [i.e both "hydrodynamic" and 
"electric" singularities], those outside K are explicitly determined by the singularities of 
x(C) = F(f(()). Therefore the derivative 

|f(/( C )) 

is meromorphic in the entire C plane, and hence it is rational with number and order of 
singular points known beforehand. This allows one to write down its explicit expression 
up to a number of indetermined coefficients. Then the conformal mapping is expressed as 

/(C) = ^~ 1 (x(0), (4.10) 

and it remains to write and solve a set of equations for location and strength of singular- 
ities. It is the case considered prevoiusly. 



4.2 Unidirectional external field 

Suppose now that 

G = H(x)=H(hz + z)), H'(x)>0. (4.11) 
It means, that the external "force" has only x-component that is independent of y. Then 

©(0 = ^(^/(0 + /'(I/O)), (4-12) 

/(C) + /* (1/C) = 2#- 1 (0(C)), C G OK. (4.13) 

The functions /(C) and /*(1/C) are analytic respectively inside and outside the unit circle. 
It is the Riemann-Hilbert problem that is solved using the Cauchy-type integral (cf.[7],[5]): 

m = 1 I g-WDj, _ i / g-(eM)^ (414) 



™ J d K U-C 27T2 J 8K U 

Example Let H(x) = x 2 , and the flow is generated by a dipole at a location 
z = x > 0. Then 

W(z) ^— , z->x . (4.15) 

z - x 

We assume that x = /(0). Then 6(C) should have poles at C = and C = oo, and, 
therefore, 

©(C)=«(c+0 +P- (E.Q) 
19 



Introducing these expressions into Eq.(4.14), we get 

= 1 / ^(u + u-p + ff iu _ 1 r V a(u + u-i) +l 3 du 
™ J 9 K U -C 2ttz J dK u 

Then the dipole location is given by the expression: 

1 / Ja(u + u~ l ) + 3 , 
x = —(p ^ >—-du, (4.18) 

while for its strength \i we find 



" -a- f(0) = ^ = l/ ^ a{u + U ~ 1)+(3 du. (4.19) 



/'(0) a 7ri J 9X u 

Equations (4.18) and (4.19) serve to find a and 3 for given x and \i. They can be reduced 
to equations 

x = — / A/2acos(/9 + /Sdv?; ~ = ~ V 2« cos <p + 3 cos </?d</?. (4.20) 

Relations Eq.(4.20) are shown in Fig. 5, a. Using these relations, we can construct explicitly 
the equilibrium domains predicted by conformal mapping Eq.(4.18). Some results are 
presented in Fig.5,b. 

4.3 Axially-symmetric external field 

We assume now that the external potential has the form 

G = H(x 2 + y 2 ) = H(zz), H'(r) > 0, r ^ in D, r = (zz) 1/2 . (4.21) 

that corresponds to an a radially-symmetric external field with the symmetry axis outside 
the equilibrium domain D. Equation (4.4) implies 

/(C)/*(l/C) = H-\e(()), (4.22) 

or 

In /(C) + In /* (1/C) = In ^(6(0). (4.23) 

As by assumption f(() ^ 0; £ e D, the logarithms in the l.h.s. of this equation are 
analytic functions respectively in the unit disk and outside it, and therefore we once more 
have the Riemann-Hilbert problem. Its solution is 

- exp (ss £ -^*' - « £ — ■ (424) 

These expressions allow us to restore the shape of the equilibrium domain provided the 
expression for 6(C) can be guessed using properties of the hydrodynamic singularities. 
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Figure 5: a: Relation between geometric parameters of the equilibrium domain 
and relative strength of the dipole; b: Shape of equilibrium domains for B = 
2.00; 2.02; 2.06; 2.12; 2.20 (curves 1-5 respectively). 
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Example \.2. Let H(x 2 + y 2 ) = r 2 , and the flow is generated by a dipole at a location 
z = r . Then 

W(z) ~ ^ * - r . (4.25) 

We assume that ro = /(0). Then repeating argument of the previous subsection, we find 
the same expression (4.16) for 0(C), and keeping in mind that in our case H~ X {X) = X, 
we have, upon introducing this expression into Eq.(4.24), 

/(CHexpf-L/ i / 

\2m J dK u-C 4rn J 9K u 

Characteristic shapes of the equilibrium domains predicted by the mapping Eq.(4.26)are 
shown in Fig. 6. 

4.4 External field depending on a harmonic function 

Let the external potential depend on a function harmonic up to specified logarithmic 
singularities, 

G = H{T{x,y)), 

M 

at = QmK* -x' m ,y- y'J- (4.27) 

m=l 

Then T(x,y) is the real part of an analytic function S(z) having specified logarithmic 
singularities, and 

G = H(±(Z(z) + Z*(z))). (4.28) 

Let function H^ 1 be rational and all "hydrodynamic singularities" correspond to multi- 
poles (there is no logarithmic singularities corresponding to sources). Then from Eq.(4.28) 

S(/(C)) + S*(/*(1/C)) = 2^(0(0), (edK, (4.29) 

or, denoting 

z(0 = s(/(0), 

Z(()) + Z*(l/C)) = H- 1 (Q(()), (edK. (4.30) 

It is essentially the same equation as Eq.(4.13), and it can be solved using the same 
technique. Then the conformal mapping 

/(C) = 5- 1 (Z(C)). (4.31) 
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Figure 6: Equilibrium domains for flow driven by a dipole in quadratic axisymmetric 
external potential field for B = 2.0; 2.0212.061; 2.121; 2.201 (curves 1-5 respectively), (a); 
blow-up of Fig. 6, a, (b). r,o 



Example 4-3. Let 

T = x 2 -y 2 ; H(T) = v ^f; E(z) = \z 2 - (4.32) 

and let the flow is generated by a single dipole of the strength /i at z = a > 0, 

W(z) — , z -> a. (4.33) 

z — a 

Then 6(C) is expressed by Eq.(4.16), and 



^ />x 1 / Jaiu + u- 1 ) + , 1 / Vafu + u- 1 ) + # , 

/(C) = v^Z. (4.35) 

Then for a and /3 we have equations 



a = /(0) = y ^ jT V 2 cos<p + P/ad<p; (4.36) 

y^y = a; /'(0) = £ = (2/Z(0)) 1 /^'(0) = 
AJa r 

/ y 2 cosy? + f3 / a cos (pdip . (4.37) 

an J 

Shapes of the equilibrium domains for flow driven by a dipole at z — 1 in the external 
field corresponding to Eq.(4.32) are shown in Fig. 7. 

4.5 Non-planar Hele-Shaw cell 

Consider a non-planar Hele-Shaw cell in constant gravity field. Let (x, y) be coordinates 
in the horizontal plane, and h(x, y) is elevation of a cell point over the horizontal plane. 
Then assuming h = h(x) it is possible to introduce the conformal coordinate 

z = s(x)+iy, s(x)= f v 7 ! + h'(t) 2 dt, (4.38) 

Jo 

so that the problem reduces to that of planar Hele-Shaw cell with effective potential 

H{z) = h(x) = h(s-\Rez)). (4.39) 

Similarly, if 

h(x,y) = K(r); r = ^x 2 + y 2 , (4.40) 
then the conformal coordinate is 

z = e llp R(r); R(r) = exp f Jl + K'(p) 2 ^, (4.41) 

Ji P 

and the effective potential is 

H(z,z) = h(RT 1 (\z\)). (4.42) 
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Figure 7: Equilibrium domains for flow driven by a dipole in external potential field of the 
form Eq.(4.32) for B = 2.0001; 2.0201; 2.0601; 2.1201; 2.2001 (curves 1-5 respectively). 
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5 Possible applications 



In this Section we show that the mathematical model considered above can be used to 
model electroosmotically-driven flow in a thin gap between to infinite parallel walls pro- 
vided that the gap is filled with two immiscible fluids having equal electric conductivities, 
viscosity of the fluid in the exterior of the domain D(t) is negligible, and electroosmotic 
coefficients of the two fluids are different. Therefore, there exists at least one non-trivial 
physical situation corresponding to the mathematical problem considered in this paper. 

5.1 Electrokinetic Effect: Physics, Available data 

Electrokinetic effect consists in generation of electric current by fluid flow through porous 
media or thin gaps between solid walls, and in the reverse effect of inducing flow by 
application of electric field. It is the last case, usually referred as electroosmosys that 
serves as primary motivation of presented theory. Electrokinetic phenomena are caused 
by difference in mobility of ions, some of which are fixed at the surface of the solid skeleton 
(matrix) of the porous medium, or the solid walls, while dissolved counterions can move 
with the fluid within the gap or porespace, or force it to move, if an electric field is applied. 
Macroscopically, the flow and electric current are governed by the equations 

u=--(Vp- W), (5.1) 
V 

I = -S(W-CVp), C=l/i. (5.2) 

Here, k is the medium permeability, r\ is the fluid viscosity, S is the fluid electric conduc- 
tivity C is the electrokinetic coupling coefficient (see [9], [2], [6] for details). 

Both "streaming potentials", i.e. electric fields generated by fluid flow, and "electroos- 
motic flow" , the flow driven by electric potential differential, have important applications. 
Electroosmotic flow is used in soil remediation and prevention of moisture penetration in 
underground structures. Recently, electroosmotic flow is also actively studied as an ele- 
ment of microfluidic devices, when flow in narrow gaps or channels is driven by electric 
potential [12]. Presumably, it is this class of flows, to which the presented above theory 
can find some applications. 

Namely, we consider flow driven both by pressure gradient and external electric field 
in a narrow plane gap between two solid non-conducting walls. We assume, that due to 
significant fluid conductivity the flow effect on electric current is negligible. In this case, 
Eqs.(5.1) and (5.2) become 

u = --(Vp- W), (5.3) 
rj 

I = -SVip. (5.4) 

here, u(x, y) and I(x, y) are averaged over the gap thickness flow velocity and electric 
current; they satisfy the continuity equations (conservation laws) 

Vu = q u (x, y); V-I = Q T (x, y); (5.5) 
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k = h 2 /12; pressure and the electric potential are functions of the in-plane coordinates 
(x,y). 

Now we assume that the gap is filled by two fluids, one of them, within time-dependent 
plane domain D(t), is characterized by the viscosity rj, conductivity S and electrokinetic 
coupling coefficients C and £; another, outside D(t), is filled by another fluid with viscosity 
rji and conductivity Si, and electrokinetic coupling coefficients C\ and £1 = \fC\. Then 
at the boundary T(t) = dD(t) we have 

p+=p- ] ^+=-0-; u+ = u~- I+ = I~. (5.6) 

For given densities of the volume (q u ) and electric (q T ) sources Eqs.(5.3)-(5.6) define a 
free boundary problem of coupled pressure/electroosmotically driven flow in the gap. 

Now we consider a particular case when both fluids have the same conductivity, Si = S. 
Then if)(x,y) satisfies the equation 

A#r,y) = | (5.7) 

in the entire plane. 

It can be considered as known. Let now the viscosity of the external fluid outside D(t) 
be negligible. Then, assuming that there is no net flux to infinity, we have 

p + (x,y) =^(x,y,t) + amst 1 , (x, y) e Z \ D(t). (5.8) 

Now we define "effective pressure" as 

P(x, y) = p(x, y) - ^{x, y) - amst x . (5.9) 

Then we have 

V- u = 0, u=--(VP-(f-6)V^); xeD(t); P(x) = 0, x e dD(t) (5.10) 
V 

It is, up to notations, the problem considered in this paper. Above examples show that 
external electric field can be used to confine the flow to a finite domain D. 

References 

[1] Varchenko, A. N. and Etingof, P. I. Why the Boundary of a Round Drop Becomes 
a Curve of Order Four. Providence, RI: Amer. Math. Soc, 1992. 

[2] Dukhin, S.S. and Deryagin B.V., Surface and Colloid Sci., Vol. 7: Electrokinetic 
phenomena, E.Matijevich, ed.; NY, 1974, J.Wiley. 

[3] Entov V.M., Etingof P.I., Kleinbock D.Ya. Hele-Shaw flows with free boundaries 
produced by multipoles, Euro Jl. Appl. Math., 4(2), 97-120, 1993 



27 



[4] Entov V.M., Etingof P.I., Kleinbock D.Ya. On nonlinear interface dynamics in Hele- 
Shaw flows, Euro Jl Appl Math, 6(5), 5, 399-420, 1995 

[5] Gakhov F.D. Boundary Value Problems, Dover, 1990, 581pp 

[6] Marino, S., Coelho, D., Bekri,S., and Adler, P.M., Electroosmotic phenomena in 
fractures, J. of Colloid and Interface Science, 223(2), 292 - 304, 2000 

[7] Muskhelishvili, N.I. Singular Integral Equations, Dover, 1992, 447pp. 

[8] Ockendon, J.R., and Howison, S.D. Kochina and Hele-Shaw in modern mathematics, 
natural science and industry, Journal of Applied Mathematics and Mechanics, 66 
(3), 505-512, 2002 

[9] Overbeek, J.Th.G., Electrochemistry of the double layer, in: Colloid Science, edited 
by H.R.Kruyt, pp.115-193, Elsevier, New York, 1952 

[10] Richardson, S. Hele-Shaw flows with a free boundary produced by the injection of 
fluid into a narrow channel, J. Fluid Mech., 56(4), 609-618, 1972 

[11] Richardson, S. Hele-Shaw flows with time-dependent free boundaries involving a 
multiply connected fluid region, Euro Jl Appl Math, 12, 571-599, 2001 

[12] Wong P.K., Wang T.-H., Deval J.H., Ho C.-M., Electrokinetics in micro devices 
for biotechnology applications, IEEE/ASME Transactions on mechanotronics, 9(2), 
366-376, 2004 



28 



